#ifndef TLS_GEOMETRY_H
#define TLS_GEOMETRY_H
#include "vector.h"
#include "transform.h"
#include <algorithm>
#include <cstdio>

namespace tls{

  struct ray{
    vec3 o,d;
    real t;

    ray(const vec3& ori,const vec3& dir,const real& s=inf):o(ori),d(dir),t(s){}
    vec3 end()const {return o+d*t;}
  };

  struct aabb{
    vec3 l,u;

    explicit aabb(const vec3& a=vec3(inf,inf,inf),const vec3& b=vec3(-inf,-inf,-inf))
      :l(a),u(b){}

    aabb make_union(const vec3& p)const{
      aabb res(*this);
      if(p.x > res.u.x) res.u.x=p.x;
      if(p.y > res.u.y) res.u.y=p.y;
      if(p.z > res.u.z) res.u.z=p.z;

      if(p.x < res.l.x) res.l.x=p.x;
      if(p.y < res.l.y) res.l.y=p.y;
      if(p.z < res.l.z) res.l.z=p.z;
      return res;
    }
    aabb make_union(const aabb& box)const{
      aabb res;
      res.u.x=cgmax(u.x,box.u.x);
      res.u.y=cgmax(u.y,box.u.y);
      res.u.z=cgmax(u.z,box.u.z);

      res.l.x=cgmin(l.x,box.l.x);
      res.l.y=cgmin(l.y,box.l.y);
      res.l.z=cgmin(l.z,box.l.z);
      return res;
    }

    vec3 center()const{return (l+u)*0.5;}

    bool in(const vec3& p)const{
      return p.x >= l.x && p.x <= u.x &&
	p.y >= l.y && p.y <= u.y &&
	p.z >= l.z && p.z <= u.z;
    }
    bool intersect(ray& r)const{
      real maxt=r.t, mint=0;
      for(int i=0;i<3;++i){
	real t0=(l[i]-r.o[i]) / r.d[i];
	real t1=(u[i]-r.o[i]) / r.d[i];
	if(t0 > t1) swap(t0,t1);
	if(t1 < maxt) maxt=t1;
	if(t0 > mint) mint=t0;
      }
      if(mint <= maxt){
	r.t=mint;
	return true;
      }
      return false;
    }

  };

  struct geom_tri{
    const vec4 *p0,*p1,*p2;
    aabb box;

    geom_tri(const vec4* a,const vec4* b,const vec4* c):
      p0(a),p1(b),p2(c),box(){
      for(int i=0;i<3;++i){
	box=box.make_union((*this)[i].v());
      }
    }


    const vec4& operator[](int i)const{
      if(i==0) return *p0; else if(i==1) return *p1; else return *p2;
    }

    vec3 normal()const{
      vec3 u=(*p1 - *p0).v();
      vec3 v=(*p2 - *p0).v();
      return cross(u,v).normalize();
    }

    bool intersect(ray& r)const{
      mat3x3 m(-r.d,(*p1 - *p0).v(),(*p2 - *p0).v());
      vec3 b(r.o - p0->v());
      mat3x3 inv=inverse(m);
      vec3 x=inv * b;

      real t=x[0],u=x[1],v=x[2];
      if(t>0 && t<r.t &&
	 u>=0 && u<=1 &&
	 v>=0 && v<=1 &&
	 u+v<=1){
	r.t=t;
	return true;
      }
      return false;
    }
  };

  struct mesh{
    struct tuple{
      int a,b,c;
    };

    int n_verts,n_faces;
    vec4* verts;
    tuple* faces;
    ~mesh(){
      delete[] verts;
      delete[] faces;
    }
  };

  void* create_mesh(std::FILE* ply);
}//namespace tls

#endif

